A first-principles approach to electrical transport in atomic-scale nanostructures 
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We present a first-principles numerical implementation of Landauer formalism for electrical trans- 
port in nanostructures characterized down to the atomic level. The novelty and interest of our 
method lies essentially on two facts. First of all, it makes use of the versatile Gaussian98 code, 
which is widely used within the quantum chemistry community. Secondly, it incorporates the semi- 
infinite electrodes in a very generic and efficient way by means of Bethe lattices. We name this 
method the Gaussian Embedded Cluster Method (GECM). In order to make contact with other 
proposed implementations, we illustrate our technique by calculating the conductance in some well- 
studied systems such as metallic (Al and Au) nanocontacts and C-atom chains connected to metallic 
(Al and Au) electrodes. In the case of Al nanocontacts the conductance turns out to be quite de- 
pendent on the detailed atomic arrangement. On the contrary, the conductance in Au nanocontacts 
presents quite universal features. In the case of C chains, where the self-consistency guarantees the 
local charge transfer and the correct alignment of the molecular and electrode levels, we find that 
the conductance oscillates with the number of atoms in the chain regardless of the type of electrode. 
However, for short chains and Al electrodes the even-odd periodicity is reversed at equilibrium bond 
distances. 



I. INTRODUCTION 

Molecular- and atomic-scale electronic devices are at- 
tracting an ever increasing interest due to the impact 
they are expected to make in the world of Nanotechnol- 
ogy. The number of experimental and theoretical works 
in this particular area of research, generically known as 
molecular electronics m, is growing exponentially. The 
design of devices at the molecular and even atomic scale 
poses new challenges which require new theoretical and 
experimental techniques to be developed. Scanning Tun- 
neling Microscopy (STM) is probably the pioneer of the 
experimental techniques in this research area. It can be 
used not only in the tunneling regime to image adsor- 
bates |^||, but also in the contact regime to build few- 
atom nanoscopic contacts M. STM can also be used to 
investigate the electrical properties of nanotubes B and 
DNA molecules || with one or both of their ends at- 
tached to a suitable electrode. In addition to STM, me- 
chanically controllable break junctions have also revealed 
themselves as powerful tools to study electrical transport 
in metallic nanobridges M or individual molecules |]||. 

The basics to calculate the zero-bias conductance G of 
a nanoscale contact had been established by Landauer 's 
in his pioneering work |lC| ] long before these systems were 
commonplace. In Landauer's formalism G is simply given 
by the quantum mechanical transmission of the electrons 
around the Fermi energy fTi]] . The value of this trans- 
mission is essentially determined by the region where the 
number of channels available for conduction is the small- 
est. In molecular or atomic-scale nanocontacts the re- 
gion of relevance is the molecule and/or the few atoms 
forming the nanoscopic bridge between electrodes. The 



transmission is thus strongly dependent on the partic- 
ular molecule, the detailed atomic arrangement of the 
electrodes in the contact region, and the chemical nature 
of them. Knowing the atomic arrangement of the elec- 
trodes or the way the molecule binds to the electrodes is, 
however, a major problem in itself. Furthermore, even if 
these important details were known, implementing Lan- 
dauer's formalism still requires to know the electronic 
structure and this is a formidable task as well. 

Calculations based on tight-binding or semi-empirical 
models |l2| jlj] have been, and still are, very popular 
since these models capture the atomic-scale character in 
some detail and are easy to implement. However, they 
do not allow for structural relaxations to be performed. 
Most importantly, these simple models, in general, do 
not yield correct values for the local electronic charges. 
In other words, the chemical potential is not uniform 
across the entire system in equilibrium. While imposing 
local charge neutrality is a straightfoward improvement 
on these models for metallic nanoconstrictions |l3| . |l4| . 
there does not exist any simple modification in the case 
of more complex systems like metal-molecule-metal het- 
erostructures p5[ , A way around this problem is to per- 
form self-consistent first-principles calculations which, at 
least at a mean-field level, guarantee the uniformity of the 
chemical potential. However, most numerical implemen- 
tations commonly used to carry out ab-initio electronic 
calculations are either restricted to finite systems, such 
as the Gaussian code or require the infinite system 
to be periodic such as the SIESTA code |jl7). None of 
these methods are suitable to address the systems here 
studied which are both infinite and non-periodic. Fi- 
nally, a perhaps more serious difficulty is the intrinsic 
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non-equilibrium character of electrical transport. 

In recent years several proposals have appeared to 
tackle this problem Jl8|-p3t. Most of all are based upon 
density functional (DF) theory. In addition to the well- 
known virtues of the DF theory, it presents the addi- 
tional advantange that Landauer's theoretical framework 
does not need to be modified since DF is still a single- 
particle description of the many-body problem. In the pi- 
oneering works of Lang and co-workers, Tsukada and co- 
workers, and Guo and co-workers |l8|,[l9],|4|-[2(| the elec- 
trodes were described within the jellium approximation. 
Jellium models are still used today ]2l|j27| , p8[ ; for they 
are convenient in one way: They provide featureless con- 
tacts which represent generic situations. As mentioned 
below this is one major feature of our model too. How- 
ever the jellium model presents serious drawbacks. How 
can the differences observed experimentally between, for 
instance, Al and Au electrodes, be taken into account by 
means of a jellium?. Furthermore, this approach is not 
satisfactory when one is trying to describe, for instance, 
STM experiments where the detailed atomic structure of 
the tip determines, to a large extent, whether or not the 
STM can resolve the topography or molecular structure 
of the adsorbate. Recent approaches, which essentially 
differ only on the numerical implementation, intend to 
incorporate the atomic structure of the electrodes in the 
DF calculation [§(|^]^|^|(J . It is pertinent noting 
here that in most of these studies a periodic structure 
beyond a given point within the leads is assumed. Effi- 
cient transfer matrix techniques |3l],[32| make this reason- 
able assumption tractable, but it forces one to consider a 
very specific type of lead (typically a finite-section wire 
p0| , p2 29 3(J or an infinite surface |2^]). Using a jellium 
model for the electrodes is harmless but it lacks the min- 
imum atomic detail which is crucial to describe, e.g., the 
contact with molecules. However, employing well-defined 
specific electrodes is not desirable either since their own 
electronic structure can interfere with the interpretation 
of the results. For instance, the appearance of gaps in 
the conductance close to the Fermi energy for perfectly 
conducting systems such as Au chains can only be con- 
sidered an artifact due to the unphysical electrode model 
|^2|,^0|. Furthermore, actual nanocontacts are not ex- 
pected to have high symmetry. As explained below these 
difficulties are circumvented in our method. 

Recently |33| , we have presented an alternative to the 
ab-initio methods above mentioned. Close in spirit to 
those presented in Refs. (2Cjj2^], the main differences and 
advantages with respect to them are the use of the stan- 
dard Gaussian98 code to carry out the DF calculation of 
the relevant transport region and the description of the 
electrodes bulk by means of appropriate Bethe lattices 
p4] , |35| . The Gaussian98 code provides a versatile method 
to perfom first-principles calculations of clusters, incor- 
porating the major advancements in the field in terms 
of functionals, basis sets, pseudopotentials, etc.. On the 



other hand, the Bethe lattices are two- fold convenient: 
(i) They reproduce the essential features of the bulk den- 
sity of states and (ii) their directional self-energies can 
be easily calculated (see Appendix). In Ref. |33| these 
ideas were applied to investigate electrical transport of a 
Ceo molecule in between Al electrodes. Here we describe 
in detail an improvement to our previous approach |35| ] 
which can be summarized in that we now incorporate 
self-consistently the semi-infinite electrodes into the con- 
ductance calculation within the Gaussian98 code. This 
requires working with Green's functions from the very 
start. The method has some resemblance with the cluster 
Bethe lattice method developed to investigate the phys- 
ical properties (electronic structure, phonons, etc.) of 
disordered systems jjji). We illustrate the possibilities 
of our method by investigating electronic transport in 
Al and Au nanocontacts. These systems have been the 
subject of extensive studies in the past by means of tight- 
binding methods [Ojljl and, more recently, by means of 
ab-initio methods ]21 ^2|,[30| . We show that, even in these 
systems where charge transfer is apparently unessential, 
tight-binding methods may fail to provide a correct quan- 
titative picture. It turns out that the conductance in 
the case of Al nanocontacts is strongly dependent on the 
detailed atomic structure. In the case of Au nanocon- 
tacts, on the contrary, the results are more universal as 
confirmed by experiments. Next, we choose a system 
with somewhat appealing features: Carbon- atom chains. 
As shown in Ref. pij these chains, when contacted by 
Al electrodes, exhibit a conductance that oscillates with 
the number of atoms in the chain. Here we address this 
problem by taking proper account of the binding to the 
electrodes and investigate how the results depend on the 
type of electrode (Al or Au). 

The rest of the paper is organized as follows. In Sec- 
tion II we discuss the main characteristics of our method. 
Section III is devoted to an extensive discussion of the re- 
sults. Finally, we end the paper by summarising the main 
features of the method and the most remarkable results 
(Section IV). 



II. THE GAUSSIAN EMBEDDED CLUSTER 
METHOD (GECM) 

In previous work [p3[ we have presented a method to 
study transport in atomic-scale and molecular devices 
which is based on standard quantum chemistry calcula- 
tions with the Gaussian98 code |l6| . This scheme, which 
has been recently adopted by other groups (see, e.g., Refs. 
p6| ) , is taken here a step further. A DF calculation of the 
region that includes the molecule or set of atoms forming 
the contact between electrodes and a significant part of 
the electrodes is performed (see Fig. |IJ). As far as trans- 
port is concerned, the hamiltonian (or Fock matrix F) of 
this central cluster or supermolecule contains the relevant 
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information since it embraces the region with the smallest 
number of channels for conduction. However, according 
to the usual theoretical transport schemes its asso- 
ciated Green's functions are unsuitable for the evaluation 
of the current (note that they simply have poles). The 
retarded(advanced) Green's functions associated with F 
needs to be extended to include the rest of the semi- 
infinite electrodes: 

G r{a) {e) = (elS — F ± iS)- 1 -> [elS — F — E r W(e)] _1 . 

(1) 



In this expression 



,r(a) . 



(2) 



where Sfl(Si) denotes a self-energy matrix that accounts 
for the part of the right (left) semi-infinite electrode which 
has not been included in the DF calculation. S is the 
overlap matrix and / is the unity matrix. The added 
self-energy matrices can only be explicitly calculated in 
ideal situations, which, in principle, limits the desired ap- 
plicability of the above procedure. For instance, in Refs. 
p2| , [30f the authors consider finite-section wires as elec- 
trodes. As a result of this choice gaps appear in the con- 
ductance of otherwise perfectly conducting central clus- 
ters. 
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FIG. 1. Schematic view of a cluster where phantom atoms 
from the Bethe lattices are shown. 

In order to overcome this type of problem, we choose 
to describe the bulk electrode with a Bethe lattice tight- 
binding model with the appropriate coordination num- 
bers and parameters (see Appendix). The advantage 
of choosing a Bethe lattice resides in that it reproduces 
fairly well the bulk density of states of any metallic elec- 
trode, avoiding this way the appearance of spurious re- 
sults. In addition to this the self-energy matrices that ap- 
pear in Eq. [j] can be calculated iteratively in a simple way 
(see Appendix for more details). For each atom of the 
outer planes of the cluster, we choose to add a branch of 



the Cayley tree in the direction of any missing bulk atom 
(including those missing in the same plane). In Fig. |l| 
the directions in which branches are added are indicated 
by smaller atoms which represent the first atom of the 
branch in that direction. Assuming that the most impor- 
tant structural details of the electrode are included in the 
central cluster, the Bethe lattices should have no other 
relevance than that of introducing a featureless reservoir. 

In our present approach the self-consistent process does 
not stop once the finite central cluster has been solved. 
Instead, we reformulate the Gaussian98 code to proceed 
with the self-consistency of the now infinite system. More 
specifically, once self-consistency for the finite cluster has 
been almost attained, we calculate the Green's function 
as explained above. Then, the density matrix is obtained 
from the Green's function according to 



Im 



G r (e) 



dr. 



(3) 



where 6f is the Fermi level fixed by the condition of over- 
all charge neutrality in the cluster. The integral in Eq. 
(^) is calculated along a contour in the complex plane 
as explained in Refs. [^,^2 23 with an efficient auto- 
matic numerical integration scheme of Perez- Jorda et al. 
p7{ . The density matrix is now used to recalculate the 
matrix elements of the Fock's operator and the process 
is repeated until self-consistency is achieved. We note 
that in this method the standard eigenvalue problem, 
inherent to the Gaussian98 code, is replaced by the cal- 
culation of Green's functions. In the end the Green's 
functions describe an infinite system in a more consis- 
tent way than in the method discussed previously by us 
p3[ since it effectively removes finite-size effects in the 
self-consistency. It is interesting to note that the appli- 
cability of this approach, which we hereafter name the 
Gaussian Embedded Cluster Method (GECM), goes well 
beyond the present study. In fact it could be a power- 
ful tool whenever an infinite media has to be described 
(for instance adsorption of molecules on solid or liquid 
surfaces). 

The conductance can now be simply calculated 
through the expression [nof 



2e 2 

G= ff_Tr[T], 



(4) 



where Tr denotes the trace over all the orbitals of the 
cluster and T is the transmission matrix which, in turn, 
is given by 



T = T L G r T R G a , 



(5) 



where the matrices Tr and Tl are given by i(EJj — SJj) 
and — fif,), respectively. In order to single out the 
contribution of individual channels to the current one 
can diagonalize the transmission matrix. It turns out 
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(see below) that only few channels give a non-negligible 
contribution to the current. The symmetry of each chan- 
nel was identified by looking at its weight on the atomic 
orbitals of the central atom of the constriction or the C 
chain. 

Finally, it is worth pointing out that if the current 
in the finite-bias regime has to be obtained, one should 
simply integrate in energy the above expression with ap- 
propriate Fermi distributions functions. Note, however, 
that a new definition of the density matrix generalized 
to non-equilibrium j2^j2^,^] needs to be used in the cal- 
culation of the Green's functions. Apart from this, the 
Landauer-type expression (|j) remains valid as long as one 
does not give up the single-particle description. In this 
work we are concerned with basic and still open aspects 
of transport in the systems studied and we will focus on 
the linear regime. 



III. RESULTS 

For all the DF calculations we have used the Becke's 
three-parameter hybrid functional using the Lee, Yang 
and Parr correlation functional (B3LYP) |38| together 
with the semilocal shape consistent pseudopotential 



(SCPP) and basis sets of Christiansen et al. |39| - |4l| . We 
have selected this combination of exchange-correlation 
functional and pseudopotential for two reasons: First, 
B3LYP is one of the most accurate and certainly the 
most popular among the gradient-corrected exchange- 
correlation functionals; second, Christiansen et al. SCPP 
provide accurate results for a wide variety of atoms while 
retaining the simplicity of a minimal basis set |39] ff2[ . 
Needless to say that there is no need to restrict the cal- 
culations to a minimal basis set since the Bethe lattice 
can be constructed for any basis set. It is only a matter 
of computational convenience that we have done it so. 
Nevertheless, in some cases, we have checked that bet- 
ter basis sets and the use of other exchange-correlation 
functionals do not modify the main conclusions of our 
work. 

In all cases we have investigated the influence of the 
number of electrode atoms included in the DF calcula- 
tion on the conductance. In general, the results do not 
vary qualitatively with the cluster size, allowing us to 
extract some general conclusions. However, at a quanti- 
tative level, this is not always the case, particularly for 
Al nanocontacts. Finally we note that, although most of 
the calculations were carried out taking the bulk inter- 
atomic distances for the electrode atoms, in some cases 
we investigated the effects of ab-initio relaxations. 



A. Metallic nanocontacts 

A complete theoretical study of electrical transport in 
metallic atomic contacts requires a realistic modelization 
of the formation process of these nanocontacts. Some 
structural studies using molecular dynamics |43|| for Al 
and ab-initio relaxations for Al [Q and Na jf5| have been 
reported. This is, however, a problem beyond the scope 
of this work. Here we consider archetypical atomic struc- 
ture models that are likely to appear in the last stages 
of the formation process of atomic contacts before the 
break-up: Single-atom contacts and atomic wires. More 
specifically, our first structure consists of two opposite 
pyramids grown in the (001) direction and "glued" by a 
single atom [see Fig. Single-atom contacts have been 
studied in the past with modified tight-binding models 
|[L3|| . We find to our surprise that, even in this simple 
case, our ab-initio results are qualitatively different from 
those obtained with these models, particularly for Al. 
Our second structure is similar to the previous one, but 
with a chain of three atoms instead of a single atom [see 
Fig. [J. Finally we have studied the same chain in be- 
tween two (111) surfaces with the chain placed on top 
of a surface atom [see Fig. 0]. A similar geometry has 
been recently studied with ab-initio techniques for Au 
and Al. Our results agree with what has been reported 
for Al using a jellium model for the electrodes fjlf, but 
not entirely with what has been obtained for Au ]22|,p0[ . 
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FIG. 2. Atomic structure of the single-atom contact model 
considered in this work. The number of (001) planes increases 
by one from (a) to (c) in both electrodes increasing the size of 
the pyramids and the corresponding atomic detail in the elec- 
trode bulk. Interatomic bulk distances have been considered 
for the whole cluster. 




FIG. 3. Atomic structure of the first atom-chain constric- 
tion model considered. The number of (001) planes increases 
by one from (a) to (c) in both electrodes increasing the size 
of the pyramids and the atomic detail in the electrode bulk. 
The distance between pyramid apex atoms is 4.8A 




FIG. 4. Atomic structure of the second atom-chain con- 
striction model considered. The number of atoms in the (111) 
electrode surface increases from (a) to (c) in both electrodes 
increasing the atomic detail of the surface. The distance be- 
tween planes is 9.lA and the positions of the atoms in the 
chain have been optimized. 



1. Aluminum 

Figure || (top panel) shows the conductance vs. energy 
for a single-atom Al contact. We have considered a 3s3p 
basis set and bulk interatomic distances. Curves (a), (b), 
and (c) correspond to clusters (a), (b), and (c) in Fig. 
||, respectively. In all the clusters the contact between 
electrodes occurs through a single atom, but the number 
of (001) planes explicitly included in the DF calculation 
for each pyramid increases from (a) to (c) (remember that 
Bethe lattices are always attached to the outer planes as 
in Fig.0). 
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FIG. 5. Top: Conductance versus energy (Fermi energy set 
to zero) of the single-atom Al contact seen in Fig. ^ for the 
three cases shown there. Bottom: Individual contribution of 
the different active conduction channels for (c). The labels 
indicate the orbital nature of the channels. The primed label 
is associated with 2nd-nearest-neighbors hoppings. 

It is impossible to know the actual atomic structure of 
the metallic contact in detail unless relaxation calcula- 
tions are performed, but we do not expect the detailed 
geometry away from the neck to be important. In fact, 
as Fig. H shows, the conductance does not change sig- 
nificantly from (b) to (c), apart from minor changes in 
the fine structure. This is a clear indication that, to 
a good extent, the conductance is determined by the 
atomic structure in the narrowest region of the neck. 
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However, from our results we see that, at least, the nine 
central atoms must be explicitly considered in this ex- 
ample. This is in contrast to the conclusions drawn in 
Ref. [|l3| by Cuevas et al. using a modified tight-binding 
model. Furthermore, the value of G around the Fermi 
level is « 3 which is remarkably different from the value 
they obtained. This discrepancy is due to a combina- 
tion of facts. First, the hopping parameters that repro- 
duce bulk properties in tight-binding models are not ad- 
equate for atoms with low coordination numbers, being 
these typically smaller than the ab-initio ones. Second, 
there are non-negligible contributions from next-to-near- 
neighbor hoppings in Al. This is clearly seen in the bot- 
tom panel of Fig. || where the contributions to the total 
current of the main individual channels for the cluster 
(c) are depicted. As shown in the figure, a channel as- 
sociated with 2nd-nearest-neighbors hoppings (labelled 
a') can give a contribution of almost one conductance 
quantum Go (Go = 2e 2 /h) at the Fermi level. The main 
contribution at comes from two degenerate p x ,p y -\ike 
channels (n x ,ir y ) that account for almost two conduc- 
tance quanta 2Go- (The z-axis has been chosen along 
the main symmetry axis of the cluster). In addition there 
are two channels that have sp z character (<x) with non- 
negligible contributions at the Fermi energy that add ap- 
proximately 0.3Go and O.lGo to the total conductance, 
respectively. None of this seems consistent with the tight- 
binding results for a similar geometry |b|. 
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FIG. 6. Same as in Fig. Q, but with three Al atoms forming 
a linear chain in between the electrodes (see Fig. 

The conductance of the three-atom chain shows a dif- 
ferent behavior from the one in the previous example and 
presents the same features for the two electrode models 
considered (see top panels in Figs. ^| and 0). For the first 



electrode model we have chosen a separation of 4.8A be- 
tween apex atoms of the pyramids. The conductance 
changes appreciably from (a) to (b) in Fig. ^], but it does 
not so from (b) to (c). (This dependence is similar to 
that in the single-atom contact). In the second case we 
have chosen a separation of 9.1 A between planes and we 
have performed an ab-initio relaxation of the chain atom 
positions (only the surface layer is included in the DF 
calculation where the number of atoms increases from 7 
to 35). Here the conductance does not depend too much 
on the number of surface atoms [see Figs. [|(a), (b), and 
(c)]. In all the cases there are oscillations as a function of 
the energy which, as the bottom panels in Figs. |^ and 
show, appear mostly in the ir channels. This is reminis- 
cent of the behavior of the transmission in a Fabry-Perot 
interferometer due to scattering at the interfaces. These 
results are similar to the ones reported in Ref. pi| where 
a jellium model was used to describe the Al electrodes. 

In the three cases considered there is no trace of possi- 
ble conductance quantization. A general trend that can 
be observed in our results for the chain is that the onset 
of the transmission through the 7r channels occurs close 
to the Fermi energy This makes G strongly depen- 
dent on small variations in the positions of the atoms in 
the chain and on the atomic structure of the electrodes 
close to the chain. This might explain why the exper- 
imental conductance histograms for Al are much more 
irregular than those for Au (see below) and other metals 
However, there are many open questions regarding 
the details of the conductance steps for Al that illustrate 
the necessity of performing both relaxation and conduc- 
tance calculations at the same time [Il6f. 
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FIG. 7. Same as in Fig. but here the pyramids have been 
substituted by (111) surface planes. The number of atoms 
included in the planes are 7(a), 19(b), and 35(c) (see Fig. 0). 
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2. Gold 

We have repeated the conductance calculations for the 
same structures considered above, but now consisting of 
Au atoms (we have used a here 5d6s6p basis set). In prin- 
ciple, only the electron of the 6s orbital is expected to 
contribute to the conductance at the Fermi energy which 
should make the analysis of conductance simpler. As Fig. 
^ shows, the conductance around the Fermi energy for the 
single-atom contact varies little from (a) to (c). We note 
again that, contrary to the tight-binding predictions Jl3[ , 
the conductance at the Fermi energy surpasses Go in all 
the curves. As concerns the contribution of the individ- 
ual channels we note that the major contribution (almost 
a conductance quantum) has sp z do character (a). Two 
degenerate channels (n x ,TT y ) of p x d\ &ndp y d-i character 
(mainly d) give around 0.25Go quanta each. This should 
be expected since the d orbitals contribute significantly 
to the density of states in bulk atoms and the number of 
near neighbors (eight) of the central atom in this cluster 
is almost the bulk coordination number of an fee struc- 
ture (twelve) . Of the two channels that give a significant 
contribution at rather high energies (above 2.0 eV) one 
has sp z do symmetry (a) and the other corresponds to 
2nd-nearest-neighbors . 
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FIG. 8. Top: Conductance versus energy (Fermi energy 
set to zero) of the single-atom Au contact seen in Fig. ^| for 
the three cases shown there. Bottom: Individual contribu- 
tion of the different active conduction channels for (c). The 
labels indicate the orbital nature of them and a' indicates 
2nd-nearest-neighbors. 

It is interesting to compare the single-atom contact 
results for Al and Au. As remarked above, whereas in 
the case of Al there was a very important contribution to 



the current at the Fermi level coming from 2nd-nearest- 
neighbors hoppings, in Au this was only appreciable at 
high energies. This cannot be understood in terms of 
the respective atomic radii which are very similar (1.43 
Aand 1.45 Afor Al and Au, respectively p8|]), but rather 
it is due to the character of the wavefunctions at the 
Fermi level in each case. Namely, while in Al the density 
of states at Ep mainly comes from p orbitals that are 
rather extended, in Au the wavefunction at that energy 
mainly has an s character and is more localised. 
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FIG. 9. Same as in Fig. ^, but with three Au atoms forming 
a linear chain in between the electrodes (see Fig. ^|). 

On the other hand, the conductance around the Fermi 
energy for the three-atom Au chain shows clearly an up- 
per limit of Go (see top panels in Figs. ^ and |l0| ), and 
it does not change qualitatively with the cluster size. 
Nevertheless, the exact value at the Fermi energy is elu- 
sive, changing by as much as 20% from cluster to clus- 
ter. We have not been able to verify whether the con- 
ductance curves for larger clusters converge to a given 
one, but all the curves present a characteristic behav- 
ior: Above the Fermi energy the conductance is fairly 
constant while below oscillates and vanishes right above 
the d channel contribution. The channel decomposition 
analysis is quite simple: A single sp z cio-like channel (a) 
contributes around the Fermi energy. Nevertheless, it is 
still difficult to explain from our results the robustness of 
the quantization observed in the experiments which 
does not deviates from Go by more than a small per- 
centage over a large range of stretching force. Large- 
scale structural studies along with conductance calcula- 
tions are also desirable here in order to make a precise 
quantitative comparison with experiments. Recent ab- 
initio works [p2|j3(| for Au nanocontacts have partially 
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addressed this problem. However, as already pointed 
out, the electrode model considered there seems to in- 
troduce serious difficulties in the interpretation of their 
conductance results. 
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FIG. 10. Same as in Fig. [], but here the pyramids have 
been substituted by (111) surface planes. The number of 
atoms included in the planes are 7(a), 19(b), and 35(c) (see 
Fig. |). 



with an odd number of C atoms p4|. However, either if 
we consider an open shell electronic structure for the C 
chain, or if the edge C atoms are bonded to a metal sur- 
face, the MOs can be filled in a different way, which may 
lead to an inversion of the aforementioned trend: that is, 
chains with an odd number of C atoms having a partially 
filled 7r shell while chains with an even number of C atoms 
having this shell completely filled. The first scenario, i.e., 
an open-shell electronic structure, is actually the ground- 
state for isolated linear C chains, while the second one 
explains the findings of Lang and Avouris for the conduc- 
tance of C chains attached to Al electrodes p4| , p7| |. In 
this second situation the details of the bonding region be- 
tween the C and metal atoms is of primary importance. 
This involves knowing the geometry of the contact, the 
C-metal surface distance, and the electronic structure of 
the metal under interest. This, however, cannot be easily 
deduced by using a jellium model to describe the leads. 



a 




B. Carbon chains 

The conductance of C-atom chains attached to Al elec- 
trodes has been calculated from first principles in previ- 
ous works by Lang and Avouris p^ , p7| . In their cal- 
culations semi-infinite jellium models were used to de- 
scribe the metal electrodes and a pseudopotential for the 
C cores. The self-consistent density functional procedure 
they used (see Refs. |^ , ^7| for details) predicted an os- 
cillatory behavior in the conductance at the Fermi level 
with maxima (minima) in those chains that had an odd 
(even) number of electrons. This oscillatory behavior dif- 
fers from that corresponding to a closed shell electronic 
structure and an sp z hybridization of linear C chains. In 
this situation each C atom added to the chain provides 
two sp orbitals that contribute to the a molecular orbitals 
(MOs) and two p orbitals that contribute to the corre- 
sponding 7r MOs. With a closed-shell electronic structure 
and an even number of C atoms there is a partial filled 
7r shell, and the chain could be considered to be intrin- 
sically a conductor. On the other hand, if the chain has 
an odd number of C atoms, the 7r and a shells are com- 
pletely filled and we would have a semiconductor. Thus, 
one would expect that chains with an even number of 
C atoms would provide higher conductances than chains 



FIG. 11. Atomic structure of a C3 chain attached to hollow 
sites in between (001) fee Al surfaces. 

For this reason, we have applied the method described 
in the previous sections to these kind of systems. More 
precisely, we have calculated the conductance for C 
chains starting from 3 C atoms attached to Al and Au 
(001) fee electrodes. The distance between C atoms was 
fixed to that used in Refs. |^ , ^7| , i.e., an equal spac- 
ing of 2.5 a.u. between C atoms. The description of the 
electrode surface was reduced to four metal atoms de- 
scribing a hollow site (see Fig. [n]) in the center of which 
the C chain was attached. The importance of the C-Al 
surface distance was analyzed by making two sets of cal- 
culations for each C chain. In the first set this distance 
was kept fixed at a reference value of 2A, while in the sec- 
ond one the distance between the C chain and the metal 
electrode was that providing the minimum energy of the 
corresponding cluster. The results arc commented next. 
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FIG. 12. Top: conductance versus energy (Fermi energy 
here set to zero) of C n linear chains (n = 3 — 5) attached to 
Al at equilibrium C-Al distance (lA). Bottom: the same at a 
distance C-Al of 2A. 

1. Al electrodes 

The discussion about the conductor or semiconductor 
character of C chains given in the preceding section can 
be easily extended to the conductance of these chains 
once they have been contacted to semi-infinite electrodes. 
Note that the oscillatory character of the conductance at 
the Fermi level with the number of C atoms n in the chain 
found in these systems can be understood as a conse- 
quence of the oscillatory nature of the conductance with 
respect to the energy (reflected in Fig. |l2|). The oscil- 
lations in G(E) come from two degenerate tt channels, 
thus explaining the maximum value of 2Gq (the only ex- 
ceptions being the C3 and C4 chains where additional a 
channels appear around the Fermi energy). The num- 
ber of peaks is given by the number of ir orbitals in the 
chain which, in turn, depends on the number of C atoms 
present in the chain. The positioning of the Fermi level 
near the center of these peaks (partially filled ir shell) or 
in between peaks (completely filled ir shell) would deter- 
mine the oscillatory trend with n. On the other hand, 
both the exact position of the conductance peaks with 
respect to the Fermi level and their average size are gov- 
erned by the distance between the edge of the C chain 
and the Al electrode. 




0.5 



FIG. 13. Top: conductance at the Fermi level of C„ linear 
chains (n — 3 — 7) attached to Al electrodes. Bottom: charge 
transferred to the C chain. 




Energy (eV) 

FIG. 14. Top: conductance vs. energy (Fermi energy here 
set to zero) of C n linear chains (n = 14, 15) attached to Al 
electrodes. Bottom: the same for n = 30, 31. 

After inspection of Fig. [l2| where we plot the G(E) 
for the two sets of calculations mentioned at the begin- 
ning of this section, it is evident that at a distance of 
2Abetween the edge C atom and the Al metal surface 
we end up with the same situation found by Lang and 
Avouris namely, maxima (minima) located at n 

odd (even) (see Fig. |l3|). Nevertheless, at the equilib- 
rium distance, where the C-Al surface distance is smaller 
(lA), and the transferred charge from the metal to the 



9 



C chain is larger, the situation has changed and now the 
maxima are located on chains with an even number of C 
atoms (the only exception being again the C3 chain) as 
Fig. [jl] shows. However, one would expect that the influ- 
ence of the bonding region, a local effect, would become 
less important as the size of the C chain increases. This is 
what is actually found in the conductance of larger num- 
ber of atoms, such as those presented in Fig. [l4|, where 
we plot G(E) for n=14,15 and for n=30,31 and a dis- 
tance C-Al surface of lA. In both cases the chains with n 
odd give larger values for G(ef) than those with n even. 
The same behavior is found for the same chains attached 
to Al electrodes at a distance of 2A, but has not been 
included in the figures for simplicity. It is worth noting 
that, strictly speaking at zero temperature, the peaks in 
G(E) and, thus, the oscillations in G(ep) are expected 
to survive in the thermodynamic limit n — > 00. This ef- 
fect will disappear when the temperature becomes of the 
same magnitude as the peaks width, where G{cf) will be 
the average value for the odd- and even-n chains. 









— c 5 

c 4 

-- c 3 


2 














1 

\ 








■ 

I; 








\ 1 



1 — * ■ 1 1 - - ^ 1 

-5 5 10 

Energy (bV) 



y A 



-5 5 10 

Energy (eV) 

FIG. 15. Top: conductance versus energy (Fermi energy 
here set to zero) of C„ linear chains (n = 3 — 5) attached to 
Au at equilibrium C-Au distance (1.5A). Bottom: the same 
at a distance C-Au of 2A. 



2. Au electrodes 

Another key point in the conductance of nanocontacts 
is the nature of the metal used in the electrodes which 
has been thoroughly discussed in previous sections. This 
is specially the case when the chemical bond between 
the molecule and the metal contact changes markedly as 
reflected in Fig. |lj as compared to Fig. [l^. There we 
show G(E) for the same systems previously analyzed for 
Al but with Au instead. When we move from Al to Au 



the bonding between the C chain and the metal surface 
weakens. This reflects in a sharpening of the conductance 
peaks and in the amount of charge transferred from the 
metal to the C chain (compare Figs. |l^ and [l3]). The 
equilibrium distance between the edge C atom and the 
Au surface is also typically larger in this case: 1.5A for all 
the chains. As for Al, there also appear two it channels, 
except for the C3 and C4 chains. However, the narrowing 
of the conductance peaks makes critical the alignment of 
the Fermi level. With respect to the oscillations of G(cf), 
we observe the same tendency irrespective of the C-Au 
distance, that is: maxima located at n odd. This is also 
a consequence of the weaker bond between Au and C and 
the lesser amount of charge transferred from the former 
to the latter. 




0.6 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 

2 3 4 5 6 7 8 

FIG. 16. Top: conductance at the Fermi level of C n linear 
chains (n = 3 — 7) attached to Au electrodes. Bottom: charge 
transferred to the C chain. 

IV. CONCLUDING REMARKS 

In summary, we have developed a methodology to self- 
consistently calculate ab-initio transport properties in 
atomic-scale systems based upon the Gaussian98 code. 
This opens the doors to the use of the most stan- 
dard quantum chemistry tool to the study of transport 
through molecules, an interdisciplinary subject of in- 
creasing interest. We have chosen to study two sys- 
tems that illustrate the capabilities of our approach: 
metallic constrictions of simple and noble metals and 
C chains with a reactive electrode (Al) and an "inert" 
(Au) electrode. In the first case, we have shown that 
assuming local charge neutrality, as commonly done in 
semi-empirical methods, may lead to qualitatively incor- 
rect results. On the other hand, the system metal-C 
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chain-metal illustrated how the chemistry of the con- 
tact may determine electrical transport and, therefore, 
the incorrectness of ascribing one or another behav- 
ior to a given molecule leaving apart the specific elec- 
trode/molecule chemistry. One should speak, instead, of 
electrical transport characteristics of the whole system, 
electrode- molecule-electrode. 



N 

£ t = J2£ n . (7b) 
«=i 

This set of 2N matricial equations has to be solved itera- 
tively. It is straightforward to check that, in cases of full 
symmetry, it reduces to the single equation discussed in 
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V. APPENDIX 

In this appendix we discuss how selfenergies for Bethe 
Lattices (BL) with no symmetry can be calculated. Sym- 
metry can be broken due to either the spatial atomic ar- 
rangement or the orbitals on the atoms that occupy each 
lattice site, or both. When no symmetry exists the self- 
energy in an arbitrary direction cannot be obtained by 
rotating that for a given direction as done in Ref. [pif . 
Instead, the following procedure has to be followed. The 
method is valid for any basis set or lattice. Let n be 
the N nearest-neighbor directions of the lattice we are 
interested in and V Ti the interatomic interaction matrix 
in these directions. The selfenergies associated to each 
direction have to be obtained from the following set of 
2N coupled self consistent equations, 



t T . = V T . Ei-Eo- (E f - T, n ) 



E fl =V n El - E - (E T - E Tl ) 



n 



vi 



(6a) 



(6b) 



where i = 1, N and fj = —t\. El is the energy times 
the identity matrix, Eq is a diagonal matrix containing 
the orbital levels, V Ti is the interatomic interaction in 
the Ti direction, and E^ and Ej* are the sums of the 
selfenergy matrices entering through all the Cayley tree 
branches attached to an atom and their inverses, respec- 
tively, i.e., 
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FIG. 17. Density of states for the Bethe lattices of Al (a) 
and Au (b) obtained with the parameters of Table |l| Con- 
tinuous lines correspond to results obtained with a spd basis 
while broken lines to those obtained with either a sp (Al) or 
a sd (Au) basis. The Fermi level was set at zero energy. 
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FIG. 18. Density of states for the Bethe lattices of Ti (a) 
and W (b) obtained with the parameters of Table II. Con- 
tinuous lines correspond to results obtained with a spd basis 
and the hep and bec lattices for Ti and W, respectively. In 
the case of Ti two further curves are given which correspond 
to the fee lattice with a spd basis (broken line) or a sd basis 
(chain line). The Fermi level was set at zero energy. 

The tight-binding parameters, which include only 
nearest-neighbors interactions, used in these calculations 
are given in Table 1. The Table reports data not only 
for the metals taken as electrodes in the present work 
(Al and Au) but also for two additional metals (Ti and 
W) commonly used in experiments and/or calculations. 
All were obtained through fittings to the electronic bulk 
band structures calculated by including second or even 
third nearest-neighbors interactions Jl9|. In the case of 
hexagonal-close-packed Ti we took as nearest-neighbors 
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six out of plane and six in-plane neighbors, as actual in- 
teratomic distances differ in less than a 2% |5(J. The 
densities of states on bulk atoms are shown in Fig. ( |i~7| ) 
and (|I8|). Although some of the features of the actual 
DOS are not reproduced (as it commonly occurs in the 
Bethe lattice approximation) the overall results are sat- 
isfactory. The major discrepancy is noted for Ti since in 
this metal the Fermi level in the crystalline case lies in 
a valley of the DOS (4^]. In the calculations reported in 
this work the sp basis was used for Al. We have described 
the electronic structure of Au by means of the spd basis 
(see Fig. ^ and Table |). In the case of Au and Ti, the 
reduced sd basis gives an excessively narrow conduction 
band. Finally, we note that in the case of Ti there are no 
major differences between the DOS for the fee and hep 
lattices. 
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TABLE I. Tight-binding parameters (in Rydbergs) used in the calculation of the density of states and selfenergies of the 
Bethe lattices for the metals commonly taken as electrodes, namely, aluminium, gold, titanium and tungsten (the last two not 
considered in the present work). Orbital on-site energies are represented by d while nearest-neighbor interactions are denoted by 
Vi. The electronic configuration taken in each case are: Al-a (Ss^p 1 ), Al-b (3s 2 3p 1 3d°), Au-a ^d^s 1 ), Au-b (Sd^s^p ), Ti-a 
and Ti-c (3d 2 4s 2 4p°), Ti-b (3d 2 4s 2 ) and W (5d 4 6s 2 6p°). Actual lattices are: aluminium and gold, face-centered-cubic, titanium, 
(a) hexagonal-close-packed and (b) and (c) face-centered-cubic, and wolframium, body-centered-cubic. The parameters were 
obtained by fitting the bulk electronic band structures given in 0|. The Fermi level ef corresponding to these parameters is 
also given. 
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